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Abstract. A model for beam customization with collimators and a range- 
compensating filter based on the phase-space theory for beam transport 
is presented for dose distribution calculation in treatment planning of 
radiotherapy with protons and heavier ions. Independent handling of pencil 
beams in conventional pencil-beam algorithms causes unphysical collimator- 
hcight dependence in the middle of large fields, which is resolved by the 
framework comprised of generation, transport, collimation, regeneration, range- 
compensation, and edge-sharpening processes with a matrix of pencil beams. The 
model was verified to be consistent with measurement and analytic estimation at 
a submillimeter level in penumbra of individual collimators with a combinational- 
collimated carbon-ion beam. The model computation is fast, accurate, and readily 
applicable to pencil-beam algorithms in treatment planning with capability of 
combinational collimation to make best use of the beam-customization devices. 



PACS numbers: 87.53.Mr, 87.53.Pb, 87.53. Uv 
1. Introduction 

In heavy-charged-particle radiotherapy with protons and heavier ions, conventional 
broad-beam systems deliver variety of volumetrically enlarged standard beams and 
an optimum one of them is chosen and customized to an individual treatment target 
(Kanematsu et al 2007). The beam customization is usually made by x-jaw, y-jaw, 
and multilcaf collimators (XJC, YJC, and MLC) and custom-made accessories such 
as a patient collimator (PTC) and a range-compensating filter (RCF) with facility- 
specific variations. For example, a PTC is always attached upstream of a RCF (Hong 
et al 1996) or is optionally attached downstream of a RCF when the MLC field is 
not satisfactorily precise for the target (Kanai et al 1999). Despite inferiority in dose 
conformity, adaptiveness, and cost and labor with the accessories, the broad-beam 
delivery systems have clinical advantages in robustness against organ motion and in 
practicality of quality assurance over dynamic beam-scanning systems (Lomax et al 
2001 and Jakel et al 2001), which will remain the same in the foreseeable future. 

In treatment planning of heavy-charged-particle radiotherapy, the pencil-beam 
(PB) algorithm has been commonly used for dose distribution calculations (Petti 1992, 
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Figure 1. Side views of a beam field showing (a) the conventional PB model, 
where uniformly generated pencil beams (arrows) generated at a collimator (black) 
travel with growing spreads (curves) to enter a patient (gray), and (b) definitions 
of the broad beam's global coordinate (z, x) and pencil beam's local coordinate 
(s,i) systems, angle 0, and position x- lso , where S, I, B, and P stand for source, 
isocenter, pencil beam, and particle, respectively. 



Hong et al 1996, Deasy 1998, Kanematsu et al 1998, Russel 2000, and Szymanowski 
2001). The PB algorithm handles a therapeutic broad beam as a set of independent 
pencil beams and superposes their doses to reproduce dose fluctuation from scatter 
in the presence of heterogeneity. The pencil beams generated with certain angular 
spread at the collimator plane will grow spatially as they travel downstream side by 
side as shown in figure QJa). 

Hong et al (1996) formulated the PB total spread at a point of interest (POI) for 
a customized beam, which is represented using the standard gantry coordinate system 
(IEC 2002-3) as 

^ (Z - Z col f + <T? ci +CTp t , (1) 

where the first term quadratically adds the spatial spread from the angular spread of 
the source size viewed from the collimator, cr STC /\z co i — z src \, in travel of the collimator- 
POI distance, \z— z co i|, and the second and third terms add subsequent scatters from 
a RCF and a patient, respectively. Their model is not exactly valid for a configuration 
with a PTC downstream of a RCF, where the constant spread <7 rc f from the RCF 
would cause an unphysical penumbra at the field edge even when the POI is in the 
proximity of the PTC. Kanematsu et al (1998) proposed an approximate model such 
that the scattering by the RCF is handled only in the angular spread for penumbra 
accuracy, ignoring the spatial spread in relatively short transport to the PTC, which 
may be formulated as 

(z - z col f + a 2 pt , (2) 

where 9 IC { is the angular spread of the scatter from the RCF at height z rc f and is 
related to the spread at the POI by er rc f = \z — z IC {\ 9 TC f in |T]). 

Those models were extensively examined against measurements in lateral 
penumbra of field edges (Hong et al 1996, Kohno et al 2004a, Akagi et al 2006) though 
there has been little quantitative discussion on behaviors other than penumbra mainly 
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due to dosimetric difficulties in the presence of heterogeneity (Kohno et al 2004b). 
In fact, both |T]) and ^ should be inaccurate because the collimator height z co \ 
would always affect the total spread er tot even in the middle of large fields, where the 
collimation should never be physically effective at all. The inaccuracy is originated 
in the model where the pencil beams continue growing in the transport from the 
collimator regardless of its influence. Consequently, the PB size could be much larger 
than the heterogeneity, which would invalidate the beam model where all the involved 
particles are assumed to receive the same interactions. 

In addition, those models can not handle combinational collimation with multiple 
collimators. In carbon-ion radiotherapy at National Institute of Radiological Sciences 
(Kanai et al 1999), tumors longer than the available field length are treated with two 
unidirectional beams involving large patient movement between the deliveries, where 
the two fields are gently patched with the field edges by an upstream collimator for 
robustness against the setup errors while the other outline edges are sharply formed by 
the downstream collimator. Besides, the upstream collimators should be optimized to 
minimize unwanted secondary radiations when the final collimator is not sufficiently 
thick. In the common practice, however, the combinational collimation is not fully 
utilized nor accurately handled in treatment planning. 

In electron radiotherapy, the phase-space theory has been rigorously applied to 
the pencil-beam algorithm to deal with large scatter of electrons (Storchi and Huizenga 
1985, Shiu and Hogstrom 1991, Boyd et al 2001, and Chi et al 2005). Their method, 
the pencil-beam redefinition algorithm, effectively restricts the size of pencil beams and 
improves the accuracy against the heterogeneity. The same idea should be applicable 
to heavy-charged-particle beams. 

In this work, we develop an accurate computational model for customization of 
heavy charged beams based on the phase-space theory in analogy with the pencil-beam 
redefinition algorithm for electron beams, examine the model computation against 
measurement and analytic estimation with an example of a customized carbon-ion 
beam, and discuss the usability in practical treatment planning. 

2. Materials and methods 

2.1. Physical models 

In the PB algorithm, a subset of particles that occupy small area in the position-angle 
phase space are handled altogether as a Gaussian pencil beam (Kanematsu et al 2006). 
The development of the pencil beam is described by the Fermi-Eyges theory (Eyges 
1948, Tomura et al 1998, and Hollmark et al 2004) with phase-space distribution of 
the involved N particles, 



where spatial variance t 2 , angular variance 8 2 , and angular-spatial covariance 9t are 
defined statistically with projected transverse position t and angle 9 as shown in 
figure Qlb) . These phase-space parameters at the current position s are related to the 
source size tx S rc and the source distance (s — s src ), or the focal-spot size and the focal 
distance in radiographic terminology, as 




(3) 
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where s src is the s coordinate of the source and e z • w — 1 is the scalar product of 
the basis vectors of the z and s axes that are nearly opposing in treatment systems. 

In transport of a pencil beam, the angular variance propagates to the covariance 
and the spatial variance with increases 

A9i = 6 2 As, At 2 " = (zM+ff 2 As) As (5) 

in distance step As while 2 stays constant in the transport. Interactions with the air 
are ignored because of the small density (» 1CT 3 g/cm 3 ) while the effects are cared by 
experimental determination of the source size <r src and the initial residual range /Jo- 
in matter of effective density p (Kanematsu et al 2003), the residual range R and 
the angular variance 9 2 are modified in step As as 

AR = -p As, AW = T As, (6) 

where T is the mean of projected scattering power T = AO 2 /As for multiple scattering, 
for the step. While the original Fermi-Eyges theory employed Rossi's formulation for 
multiple scattering (Eyges 1948), Highland (1975) introduced a logarithmic correction 
term for non-stochastic influence of single scattering, which was further generalized 
and extensively tested by Gottschalk et al (1993). 

A differential form of the Highland-Gottschalk formula gives the scattering power 
for a particle with charge Z e, mass A u, and kinetic energy E as 

rr- ( 1, A 2 /14.1McV \ 2 1 „ [ s As' 

where the momentum- velocity product pv = E (E + 2 Au)/{E + Au) is quantified 
with the range-energy relationship R = R(Z,A,E) « (A/Z 2 ) R(l,l, E/A) (ICRU 
1993) and the normalized radiation length I is the distance measured in the material- 
specific radiation length Xq cumulatively along the beam down to the current point. 

Use of geometric mean pv — ^/po^oPi^i of the values before (po^o) and after 
(pivi) the step for the pv in (J7]) leads to the mean scattering power T for the step As 
at a precision of 3% or better for p As < 0.7R (Gottschalk et al 1993). Steps thicker 
than that should be recursively subdivided until the absolute error becomes small 
enough. While the effective scattering point for the step is ideally As/y/3 upstream 
of the step end (Yao et al 2006), the concurrent energy loss moves the point slightly 
downstream, approximately at the center of the step (Gottschalk et al 1993). 

Loss of primary particles and yield of secondary particles in nuclear interactions 
are implicitly and approximately involved in the depth-dose curve being referenced 
in dose calculation. This framework of the phase-space theory for charged-particlc 
beams interacting with matter is also applicable, and in fact has been applied, to dose- 
distribution calculations in patient body with fine steps to deal with the heterogeneity 
(Kanematsu et al 1998, 2006, and 2008 and Akagi et al 2006). 

2.2. Beam development models 

2.2.1. Generation A matrix of pencil beams is generated on a plane at height z = z 
in the global coordinate system, where a broad beam is first effectively modified to 
form a lateral structure by either a collimator or a RCF. The pencil beams are defined 
to have traveling distance s = in the local coordinate systems shown in figure [11(b) . 
They are placed at grid points (xoj,yoi) for row i column j with spacing So, which 
coincide with the field grids (xi SO j,yi soi ) with spacing <5i so on the isocenter plane in the 
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beam's eye view, and are comprised of N particles with residual range R and slopes 
a x = dx/dz and a y = dy/dz. These parameters are handled as a function of height z 
in transport with initial values 

x{z ) i: j = x 0j , y(z )ij = y 0i , a x (z )ij = — ; — , a y (z )y = 



zo — Z SIC Zq — z. 



src 



U ( \ n N ( z o)i] ^ ( \ S XQj y 0z Zq - Zsrc fo s 

ft\Zo)ij = -K , J2 = *iso [Xisoj, 2/iso J , 7 — — ~ = ~ — — ~ ; — ) \°) 

0i s 



^iso ^-src 



where the relative number of particles N reflects the original broad beam fluence 
expected on the isocenter plane, $i so . For the pencil beams generated with uniform 
distribution in square area Sq and with the preserved focal distance, the phase-space 
parameters are initialized to 

0*(*o)« = , alc , 2 , ei(zoh = Tn~\ ~ [' = ?l ( 9 ) 

(z -Z src )^ 12|z -Z src | 12 

where the source size cr src reflects the effects of ignored scattering materials in 
the system as well as the initial beam emittance. The beam field is conveniently 
represented with matrices of the PB parameters, x, y, a x , a y , R, N, 2 , 0t, and t 2 , 
that will develop in transport as a function of height z. 

2.2.2. Transport In transport of pencil beams in the air from height zo to the next 
device at height zt, the positions and the variances, x~ij, y^, 0£y, and t 2 ^ , are modified 
by 

As 



Axij = a^ij Az, Ay^ = a^. Az, Az = z\ - z , ^- = -y a^- + a^ 2 . + 1, 

A9i lj = iPij As, A¥ij = heiij + ff 2 ^ As) As, (10) 

according to the formulation in section |2"TT1 We denote here the height immediately 
before the device as ht\, where the effects of the device have not yet been included, 
and the PB parameters are modified in the transport as t 2 (}zi)ij — t 2 (zo)ij + At 2 ij, 
for example. 

2.2.3. Regeneration Each of the transported pencil beams at height "zi immediately 
before the device contains particles whose probability density in slopes and positions 
(a x , a y , X, y) is described by phase-space distribution 



F tj (a x , Oy, x, y) = — \ {e 2 lj t 2 y - 6t t ^ 



e 2 . . t 2 . ._st?. %i 02. . t 2.._fl7?. 2 /-i-rx 

xe »3* « >3 e 13 IJ l 3 ( (ll) 

which is derived from two-dimensional extension of ([3]) with approximation 9 rs 
tan 6* (ay — a). In the absence of residual range variation, these particles are 
indistinguishable among the pencil beams and are superposed to form a broad beam 
with flucncc, mean slopes, and angular variance of the particles, 

da x / da y ^ F %j (a x , a y , x, y) = ^ -L e 2t2 *j , (12) 

-oo J-oo j • — 2 7T t ij 
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9 2 (x,y) 



®(x,y) 
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oo J —oo 
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a x 2 (x,y) +a y 2 (x,y) 



1 AT (x-x iJ ) 2 + (y- I<i j) 2 



,| " : -'--// : ^7 :>-/ 2 , ( 



t 1 
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(15) 



as a function of position (x, y) on the h\ plane, respectively. The reconstructed broad 
beam is redivided into new pencil beams at reinitialized grids on the z\ plane as 

Si _ Xij _ y u zi - Z src 

^iso ^iso^' IJisoi ^-iso — ^src 



x(zi) ij =x lj , y(zi)ij=yii, 



12 12 \Z\ — Z src 



R(zi)ij = R(zo)ij. 



(16) 



Denoting any parameter p at the height immediately before the device as p" = p(}z\), 
the statistical parameters of the regenerated pencil beams are redefined as 



N( Zl ) 



13 
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(17) 
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where the exponential, error, and hyperbolic functions will be quickly enumerated 
with the standard math library. 

2.2.4- Collimation A collimator is described as a matrix of transmission factors, 
T, with the same grids as those for the regenerated pencil beams at its downstream 
face ignoring the collimator thickness (Kanematsu et al 2006). For the pencil beams 
transported to and regenerated at z\ , the collimation modifies the relative number of 
particles immediately after the collimator at height Z\, 

NfzJi^NizJyTy, (21) 

where transmission Tij has value 1 when grid ij is in the collimator aperture or 
otherwise 0. A series of the processes of transport, regeneration, and collimation can 
be repeated for multiple collimators as long as the residual ranges have not been varied 
in the field. 

2.2.5. Range compensation A RCF is a sculptured object designed to absorb 
extra ranges of the incident particles beyond the target, which also inevitably adds 
nonuniform scattering that deteriorates particle equilibrium and consequently field 
uniformity. The RCF made of a material of effective density p and specific radiation 
length X is assumed to have a flat downstream face located at height Z\ and a 
shape described by a matrix of thicknesses, S. Pencil beam ij has a path length of 
approximately Sij in the RCF, ignoring the small beam-divergence effect. 

We first transport the pencil beams to the downstream face of the RCF at z\ 
ignoring the interactions with matter as described in section 12.2.21 which leads to 
8t(hi)ij and t 2 (*zx)ij. For the first RCF before which all the pencil beams have the 
same residual range, we regenerate them as described in section [2. 2. 31 Thus we obtain 
@ 2 (zi)ij, 9t(zi)ij, and t 2 (zi)ij. On exit from the RCF, the range Rij and the angular 
variance 9 2 ij are modified by 

ARij = -pSij, Any- 
where Z and Wiij are the charge and the mean momentum-velocity product of the 
particles in beam ij. The growths in spatial variance and covariance for the travel 
from the effective scattering point that is approximated to the midpoint, 

AOlij = AWij (0.5 Sa) , At?ij = A¥ij (0.5 S l3 ) 2 , (23) 

are correctively added to Bt(z\)ij and t 2 (z\)ij for the pencil beams exiting from the 
RCF at b zi. Note that this scattering correction along with the in-air transport is 
mathematically equivalent to a sequence of transport to the effective scattering point , 
effective-point scattering, and transport to the RCF reference face. 

2.2.6. Edge sharpening When a PTC is placed downstream of a RCF, a sharp-edged 
field must be formed by the PTC, whereas the modulated range loss and scattering 
effects originated by the RCF prevents from applying the regeneration technique. 
For an approximate solution, we partly modify the transported and collimated pencil 
beams near the collimator edge at height z\ in such a way that the spreads of the 
outgoing pencil beams at height "zi are conditionally scaled to the distances of closest 
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Figure 2. Calculated profiles in fluence $ relative to the open beam fluencc <J>o 
on the (a) collimator and (b) isocenter planes with lateral position x in units of 
local PB interval 8. The thin solid and thick gray lines are the original and ideal 
ones while the dashed and dotted lines are edge-sharpened ones with a = 1 and 
2, respectively. 



approach to the collimator edge, , with control parameter a, 



t*( b zi)« = min [% , MfzJii = , t2(bZl) % (24) 

\ Oi / \Z\ — z src \ 



a, = nun \ \ X H' ~ x n\ — rr 



" v(i'j')e(n'j'=o) 



2 / r \ 2 

01 



. 1 t vn'-yii\ - TT 

3' ¥=3 \ ^ 



i' 



(25) 



where the conditional terms are applicable when the blocked pixel i'j' is in different 
row or column. The beam spots are regulated so that the radii of a standard deviations 
will be within the aperture. 

Fi gure [2] shows an example with 7V(zi)ij — 1, T%j — H(xnj), fi[z\)%j — 5^, 
and t 2 (zi 80 )ij — t 2 (zi)ij + (5; 2 so , where the step function, H(x) = 1 for x > or 
otherwise 0, represents a half-field collimation. In this case, the edge sharpening 
successfully reduced the spread near the collimator edge though with strong distortion. 
In penumbra behavior on the isocenter plane, the edge-sharpened profiles agreed well 
with the ideal one, while the unwanted distortions mostly collapsed in the transport. 
The distortions should naturally have spatial structure as small as the spread being 
sharpened, which is normally much smaller than the subsequent spread in the patient, 
and will generally collapse in dose distributions. 



2.3. Implementation and validation 

2.3.1. Experimental apparatus We examine here how well the formulated beam- 
customization model can handle combinational collimation with one of the therapeutic 
beam lines of accelerator facility HIMAC at National Institute of Radiological Sciences 
(Kanai et al 1999). A 12 C 6+ beam of per-nucleon kinetic energy E/A — 350 MeV was 
broadened with the wobbler-scatterer system to form a 15-cm</> uniform (nominally 
±2.5%) field with a source at height 950 cm from the isocenter. The resultant residual 
range 19.6 cm in water corresponds to E/A — 327 MeV. As shown in figure^ the 
beam was customized with an XJC at 117 cm, a YJC at 96 cm, a MLC with 23 pairs 
of leaves at 69 cm, a RCF at 35 cm, and a PTC at 22 cm in height. The XJC and 
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Figure 3. Illustration of the modeled beam-customization devices in (a) side 
view and (b) beam's eye view, where the filled areas represent XJC, YJC, MLC, 
and PTC from upstream to downstream, and the hatched area represents RCF. 



YJC apertures were set to (—3.0, +3.0) cm in x and (—5.0, +3.0) cm in y, respectively, 
and the 11 lower and 12 upper leaf pairs were set to (—5.0, +5.0) cm and (—2.0, +2.0) 
cm in x, respectively. A 3-cm thick PMMA plate that only covered the left half of 
the field and an 8-cm-square aperture block were placed as the RCF and the PTC, 
respectively. 



2.3.2. Measurement Dose profiles in the air on the isocenter plane were measured 
with a 15-mm 3 2-mm0 pinpoint ionization chamber at intervals of 1 to 2 mm with 
movement precision of 0.1 mm. The profiling axes were in x at y — — 2.0 cm and +1.0 
cm and in y at x = —1.0 cm and +1.0 cm with 0.5-mm alignment precision, where the 
field edges were approximately formed by the individual collimators with and without 
the RCF. The measured doses were then divided by the corresponding open beam doses 
Do without collimation nor compensation to correct the non-uniformity of the broad 
beam. The penumbra sizes are derived from the dose-ratio profiles with dosimeter-size 
correction, tfeo^so = Vi^'la^sa ~ 1-68 2 CTj os ), where is the observed 20%-80% 

distance and <7dos = 0.5 mm is the geometrically estimated rms dosimeter size for the 
2-mm0 diameter. For dosimetric analysis, the tissue-air ratio for the 3-cm PMMA 
plate was measured with the same beam to be 0.951. 



2.3.3. Analytic estimation In the local broad-beam approximation near the 
individual collimator edges, penumbra behaviors in relative fluence with respect to 
the open beam fluence, <&/<I>0: are calculated as 

^= 2 + 2-^72^' (26) 
where the signed distance to the closest point on the collimator edge, cZ m ; n , is positive 
(negative) in (out of) the field (Hong et al 1996) and the total rms spread Ctot is 
related to the penumbra size by d2o-»80 = 1.68 otot- 

The source size is estimated to be <7 src = 2.54 cm inversely from (fTJ) with the 
measured penumbra size c?2(W80 = 4.8 mm at the upper y edge at x = +1 cm formed 
by the YJC with \z — z co i| = 96 cm and a IC { — a pt = 0. The 3-cm PMMA plate at the 
RCF adds angular spread rc f = a YC [/\z — z rc {\ = ^/(AO 2 ) = 3.3 mrad from ((6]) and 
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Figure 4. Computed fluence distribution on the isocenter plane, where the iso- 
fluence lines are 20%, 50%, 80%, 97.5%, and 102.5% and the dashed lines indicate 
the profiling locations with figure-part symbols in figure [5] 



(O with effective density p = 1.16 and specific radiation length Xq — 34.07 cm (Yao 
et al 2006). 

Equation ([I]) leads to the penumbra sizes for the two left edges in the x profiles 
and the upper edge in the y profile at x = — 1 cm, where the RCF is downstream of 
all the active collimators, ([2]) does for the lower y edge at x — — 1 cm, where the PTC 
is downstream of the RCF, and both (TTJ) and @ equivalently do for the edges in the 
right half field with no actual RCF. In dosimetric analysis, the tissue-air ratio 0.951 
for the PMMA is multiplied to the analytic fluences in the x < region. 



2.3.4- Model computation In the framework described in section 12.21 pencil 
beams were generated at the XJC, immediately collimated, transported to the YJC, 
regenerated, collimated, transported to the MLC, regenerated, collimated, transported 
to the RCF, regenerated, range-compensated, transported to PTC, collimated, edge- 
sharpened, and transported to the isocenter plane, where the grid spacing on the 
isocenter plane and the edge-sharpening parameter were chosen to be 8i so = 0.1 cm 
and a = 2 in the model computation in addition to the common parameters in the 
analytic estimation. 

The fluence distributions at the devices were calculated with ([12")) and similarly 
the in-air dose distribution on the isocenter plane was calculated with 

AT n (x-- ij ) 2 +(»-»jj) 2 

^) = E^v^ e " 2 ^ 3 ' (27) 

i j Zirt ij 

where DBBij is the tissue-air ratio amounting to 0.951 for pencil beam ij with Tij < 
or otherwise 1. 




Figure 5. Profiles in dose D relative to the open beam dose Do, along x axis 
at (a) y = —2 cm, (b) y = +1 cm, along y axis at (c) x = —2 cm, and (d) 
x = +1 cm, where the black lines, the thick gray lines, and the open circles are 
the model-computed, analytic, and measured ones, respectively. 



3. Results 

3.1. Field-formation process 

Figure [4] shows the fluence distribution on the isocenter plane in the model 
computation, which are consistent with expectations such that the field edges formed 
at the collimators become gentler with the traveling distance. The dip and the bump 
around x w are due to lateral particle disequilibrium caused by the the PMMA 
half plate at the RCF. The whole computation took only about ten seconds with 
FORTRAN interpreter in analysis package PAW by CERN on PowerPC G5 2-GHz 
processor by Apple/IBM. The short computational time is a great advantage of the 
deterministic calculation compared to Monte Carlo simulations which may need orders 
of magnitude more time (Paganetti et al 2004). 

3.2. Dose-profile analysis 

The profiles of dose ratios between the customized and open beams, D/Do, are 
plotted in figure [5j where non-uniformity of the actual broad beam should have 
been compensated. The measured profiles involve spatial uncertainties from the grid 
interpolation (0.2 mm) and the global misalignment (0.5 mm) and dose uncertainty 
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Table 1. Experimentally measured, analytically estimated, and model-computed 
penumbra sizes in the dose profiles shown in figure [5] 



Profiling 


Interested 


Effective 


20%- 


-80% penumbra 


size (mm) 


position 


edge side 


device(s) 


measured analytic 


computed 


y = —2 cm 


x left 


XJC+RCF 


6.4 


6.3 


6.0 


y = —2 cm 


x right 


XJC 


5.8 


6.0 


5.9 


y = +1 cm 


x left 


MLC+RCF 


4.6 


3.9 


4.1 


y = +1 cm 


x right 


MLC 


3.7 


3.3 


3.5 


x = — 1 cm 


y lower 


PTC+RCF 


2.3 


1.6 


2.3 


x = — 1 cm 


y upper 


YJC+RCF 


5.6 


5.2 


5.2 


x = +1 cm 


y lower 


PTC 


1.4 


1.0 


1.5 


x = +1 cm 


y upper 


YJC 


4.8 a 


4.8 a 


5.0 



a The YJC penumbras were calibrated to determine the source size in the model. 



of 0.3% from the dosimeter resolution. 

3.2.1. Penumbra Table Q] summarizes the resultant 20%-80% penumbra sizes, where 
the measured ones involve the dosimeter-size correction with the estimated uncertainty 
of 0.2 mm. The discrepancies among the measured, analytic, and computed ones 
turned out to be at a submillimeter level and are consistent with the uncertainties 
in the experimental and theoretical systems, excluding the global misalignment that 
should not be influential to the penumbra sizes. 

3.2.2. Collimator scatter In dose ratios between the customized beam and the open 
beam, the measured ones turned out to be larger than the computed counterparts by 
about 2% throughout the field. This irreproducible dose excess may have come from 
particles hard-scattered by the collimators (Kusano et al 2007a) ignored in the model. 
The collimator-scatter contribution should naturally attenuate with depth, which is 
consistent with the observation such that the excess was smaller with the PMMA half 
plate at RCF in the x < region in figure [5] 

3.2.3. Scatter modulation Scatter by the PMMA half plate caused a dip and a 
bump at x in figure [H^a) and figure [^b) , where the computed bump/dip ratios 
were (a) 1.30 and (b) 1.30 while the measured counterparts were (a) 1.28 and (b) 
1.28, respectively. Considering the dosimetric limitations, that may be an excellent 
agreement and indicates accurate evaluation of the RCF scatter in the field, which the 
conventional PB algorithms must fail to achieve. 

3.2.4- Multiplicative collimation In the model computation, the PTC aperture, 
\x\, \y\ < 4.09 cm on the isocenter plane, naturally influenced the field edges formed 
by the XJC at x — ±3.42 cm and by the YJC at y = +3.34 cm, and suppressed the 
dose tails in figure O Such effect may physically exist, but unfortunately resulted in 
larger discrepancy from the measured data for the ignored collimator scatter. 

4. Discussion 

Hong et al (1996) suggested using the height of the nearest collimator for z co \ in (TTJ) 
for combinational collimation in the PB algorithm in analogy with the broad-beam 
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algorithm. However, the nearest collimator would suddenly switch from one to the 
other in the middle of the field, where the PB algorithm would cause artifactual dose 
fluctuation regardless of heterogeneity. In contrast, the edge-sharpening process in this 
work only applies to the field-edge region and therefore there will be no unphysical 
collimator dependence in the middle of the field. 

Though we assumed here a single RCF with a flat downstream face, the 
formulation is already capable to handle multiple RCFs of any structure as well. 
In the absence of a RCF, the regeneration technique should be applied on the 
entrance to the patient to remove the unphysical collimator-height dependence. The 
collimator thickness effect can be handled in this framework with virtual multiplicative 
collimation with spacing corresponding to the thickness in a manner similar to 
and more sofisticated than the method by Kanematsu et al (2006). The present 
computational model for beam customization will provide a set of pencil beams to 
various PB algorithms to calculate dose distributions in treatment planning. 

The hard-scattered secondary particles are generally out of the scope of PB 
algorithms including this work. Practical modeling of collimator and phantom scatters 
in nuclear interactions may have yet to be studied (van Luijk et al 2001, Pedroni et 
al 2003, and Kusano et al 2007b), or could possibly be only resolved by Monte Carlo 
methods (Paganetti 2002, Kase et al 2006, and Titt et al 2008). Fortunately, the 
clinically relevant dose in ion-beam therapy is generally dominated by the primary 
particles (Matsufuji et al 2003). 

The regeneration technique here only considers particle flow in the transverse 
plane in the absence of heterogeneity. In contrast, the pencil-beam redefinition 
algorithm, which is based on the same principle, deals with electron transport in the 
presence of heterogeneity by introducing additional particle flow in the energy space 
(Shiu and Hogstrom 1991). Such extension should be also valid for heavy charged 
particles and would further improve the accuracy of dose distributions. However, 
an order of magnitude smaller scatter and an order of magnitude better spatial 
accuracy generally required for heavy-charged-particle radiotherapy could possibly 
make implementation of the energy flow less significant or less practical. 

5. Conclusions 

The phase-space theory for beam transport has been successfully applied to 
computational modeling of beam-customization devices in hcavy-charged-particlc 
radiotherapy to accurately deal with multiple scattering in the presence of multiple 
collimators and a range-compensating filter in fluence distributions at a submillimeter 
level. 

The present computational model is efficient and readily applicable to pencil- 
beam algorithms for treatment planning and will enable combinational collimation 
and compensation to make the best use of the beam-customization devices in clinical 
practice. 
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